예측모형 - 헬로 월드 에서 작업한 데이터와 모델을 가져온다.
library(tidyverse)
library(tidymodels)
library(timetk)
library(modeltime)
# 데이터 ----
full_tbl <- read_rds("data/full_tbl.rds")
# 모형 ----
wkfl_fit_lm <- read_rds("data/wkfl_fit_lm.rds")예측할 데이터와 모형개발에 활용할 데이터로 나눈 후에 모형개발에 활용할 데이터를 훈련/시험 데이터로 나눈다.
## 예측 데이터와 모형개발 데이터로 분리
forecast_tbl <- full_tbl %>%
filter(is.na(확진자))
history_tbl <- full_tbl %>%
filter(!is.na(확진자))
## 훈련/시험 데이터 분할
splits <- history_tbl %>%
time_series_split(date_var = 날짜,
assess = 30,
cumulative = TRUE)
splits<Analysis/Assess/Total>
<309/30/339>훈련/시험 데이터로 잘 나눠졌는지 시각적으로 확인한다.
splits %>%
tk_time_series_cv_plan() %>%
plot_time_series_cv_plan(.date_var = 날짜,
.value = 확진자)tidymodels 생태계의 recipe 팩키지 recipe() 함수를 사용해서 피쳐 공학(feature engineering) 작업을 수행한다. 추후, 다양한 피처를 작업할 예정이라 대략적인 틀만 잡아둔다.
recipe_spec <- recipes::recipe(확진자 ~ ., data = training(splits))
recipe_spec %>% prep() %>% juice()# A tibble: 309 x 2
날짜 확진자
<date> <int>
1 2020-01-23 0
2 2020-01-24 1
3 2020-01-25 0
4 2020-01-26 1
5 2020-01-27 1
6 2020-01-28 0
7 2020-01-29 0
8 2020-01-30 0
9 2020-01-31 7
10 2020-02-01 1
# ... with 299 more rowsparsnip 팩키지의 linear_reg() 를 통해 선형회귀(?) 모형으로 앞서 작업을 했다면, 이번에는 시계열 모형의 GLM과 같은 위치를 하고 있는 ARIMA 모형을 특정하여 예측모형을 개발한다. 이를 위해서 workflow() 를 사용하는데 기본적으로 피처 공학의 recipe_spec이 필요하고 model_spec의 모형이 필요하다.
model_spec <- arima_reg(
mode = "regression") %>%
set_engine("auto_arima")
wkfl_fit_arima <- workflow() %>%
add_recipe(recipe_spec) %>%
add_model(model_spec) %>%
fit(training(splits))
wkfl_fit_arima== Workflow [trained] ==========================================================
Preprocessor: Recipe
Model: arima_reg()
-- Preprocessor ----------------------------------------------------------------
0 Recipe Steps
-- Model -----------------------------------------------------------------------
Series: outcome
ARIMA(0,1,1)
Coefficients:
ma1
-0.2198
s.e. 0.0519
sigma^2 estimated as 2729: log likelihood=-1654.94
AIC=3313.87 AICc=3313.91 BIC=3321.33workflow 객체를 modeltime 팩키지 modeltime_table() 함수에 넣어 객체로 만든 후에 modeltime_accuracy() 함수를 사용해서 모형 성능을 파악한다.
model_tbl <- modeltime_table(
wkfl_fit_lm,
wkfl_fit_arima
) %>%
update_model_description(.model_id = 1, "선형회귀") %>%
update_model_description(.model_id = 2, "ARIMA")
model_tbl %>%
modeltime::modeltime_accuracy(testing(splits))# A tibble: 2 x 9
.model_id .model_desc .type mae mape mase smape rmse rsq
<int> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 1 선형회귀 Test 637. 77.6 6.45 128. 682. 0.825
2 2 ARIMA Test 275. 30.3 2.78 36.7 341. NA 시험데이터를 통해 모형으로 예측한 값을 시각화한다. MAE 값은 많이 개선되었으나 ARIMA 예측모형에 대한 \(R^2\) 값은 상수라서 계산이 되지 않는 아쉬움이 있다.
calibration_tbl <- model_tbl %>%
modeltime_calibrate(
new_data = testing(splits)
)
calibration_tbl %>%
modeltime_forecast(
new_data = testing(splits),
actual_data = history_tbl,
conf_interval = 0.10
) %>%
plot_modeltime_forecast(
.legend_max_width = 60,
.legend_show = TRUE,
.conf_interval_show = TRUE,
.conf_interval_alpha = 0.20,
.conf_interval_fill = "lightblue",
.title = "코로나19 확진자 1개월 예측"
)앞선 모형은 history_tbl 을 훈련/시험 데이터에 대해 적합을 시킨 것이라 … 이제 시간을 확대하여 modeltime_refit() 함수를 사용해서 모형을 전체 데이터에 대해 다시 적합시킨다.
refit_tbl <- calibration_tbl %>%
modeltime_refit(data = history_tbl)마지막으로 앞서 구축된 모형을 바탕으로 현재까지 입수된 데이터를 바탕으로 한달 후를 예측한다.
refit_tbl %>%
modeltime_forecast(
new_data = forecast_tbl,
actual_data = history_tbl,
conf_interval = 0.3
) %>%
plot_modeltime_forecast(
.legend_max_width = 25,
.conf_interval_fill = "lightblue",
.conf_interval_alpha = 0.7,
.interactive = TRUE
)마지막 단계로 데이터와 모형을 저장하여 다음 단계로 나가기 위한 작업을 수행한다.
corona_fcst_list <- list(
# 데이터 ---
data = full_tbl,
# 예측모형 ----
model = list(
wkfl_fit_lm = wkfl_fit_lm,
wkfl_fit_arima = wkfl_fit_arima
)
)
# 모형 + 데이터 저장 ----
corona_fcst_list %>%
write_rds("data/corona_fcst_list.rds")데이터 과학자 이광춘 저작
kwangchun.lee.7@gmail.com